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ABSTRACT 

An important and widely neglected aspect of the interaction between an accretion disc 
and a massive companion with a coplanar orbit is the vertical component of the tidal 
force. As shown by Lubow, the response of the disc to vertical forcing is resonant at 
certain radii, at which a localized torque is exerted, and from which a compressive 
wave (p mode) may be emitted. Although these vertical resonances are weaker than 
the corresponding Lindblad resonances, the m = 2 inner vertical resonance in a binary 
star is typically located within the tidal truncation radius of a circumstellar disc. 

In this paper I develop a general theory of vertical resonances, allowing for non¬ 
linearity of the response, and dissipation by radiative damping and turbulent viscosity. 
The problem is reduced to a universal, non-linear ordinary differential equation with 
two real parameters. Solutions of the complex non-linear Airy equation are presented 
to illustrate the non-linear saturation of the resonance and the effects of dissipation. 
It is argued that the m = 2 inner vertical resonance is unlikely to truncate the disc 
in cataclysmic variable stars, but contributes to angular momentum transport and 
produces a potentially observable non-axisymmetric structure. 

Key words: accretion, accretion discs - binaries: close - celestial mechanics hy¬ 
drodynamics - waves. 


1 INTRODUCTION 

1.1 Horizontal and vertical tidal forces 

Accretion discs in binary stars and protoplanetary systems are subject to periodic tidal forcing from the orbiting companion 
objects. Similar effects can occur in planetary ring systems, where the planetary satellites provide the perturbing forces, and 
in galactic discs. Even in the simple case of a companion with a circular orbit that is coplanar with the disc, the effects of the 
tidal force are manifold. 

Consider a system of two spherical bodies of masses Mi and M 2 in a bound orbit. Consider a third, test body orbiting 
about Mi. In a non-rotating coordinate system centred on Mi, the force per unit mass experienced by the test body is 

GMi r GM 2 (r-d) GM 2 d 

|r | 3 “ \r-d\ 3 ~ |d | 3 ’ W 

where d(t) is the position vector of the companion M 2 with respect to Mi, and the third term is the fictitious force arising 
from the acceleration of the origin of the coordinate system. Inasmuch as the second and third terms perturb the Keplerian 
motion of the test body about Mi, they represent the tidal force in this system. 

Consider a thin accretion disc around Mi, and let ( R,cf>,z ) be cylindrical polar coordinates such that 2 = 0 is the mid¬ 
plane of the disc. Suppose that the companion has a circular, coplanar orbit of radius d and angular velocity ui. Then the 
components of the tidal acceleration are 

ftR = —GM 2 [R 2 + d 2 — 2dRcos(4> — uit) + z 2 ]~ 3 ^ 2 [R — dcos(4> — uit)] — GM 2 d~ 2 cos(0 — uit), 

ftrf, = — GM 2 [R 2 + d 2 — 2dRcos((j) — cot) + z 2 ] 3 ^ 2 d sin(0 — uit) + GM 2 d~ 2 sin(<^> — uit), 

f tz = -GM 2 [R 2 + d 2 - 2dRcos((p - uit) + z 2 ]~ 3/2 z. 
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The horizontal components of the tidal acceleration are nearly independent of 2 within the disc and cause a general 
non-axisymmetric distortion of its streamlines and surface density distribution (Papaloizou & Pringle 1977). Indeed, the 
streamlines of the disc correspond closely to a family of simple periodic orbits of the restricted three-body problem found by 
Paczynski (1977). The tidal distortion is enhanced in the neighbourhood of Lindblad resonances, radii at which the forcing 
frequency resonates with the natural epicyclic oscillations of the disc. Usually at a Lindblad resonance, a non-axisymmetric 
wave will be launched and will then propagate some way radially through the disc before dissipating and transferring its 
angular momentum to the disc (Goldreich & Tremaine 1979). The resonant torque exerted between the companion and the 
disc may be calculated from a standard formula that applies under quite general circumstances in linear theory, for example 
if the response is dominated by viscosity and no wave is emitted (Meyer-Vernet & Sicardy 1987), or if the vertical structure 
of the disc is taken fully into account (Lubow & Ogilvie 1998). 

Tidal torques usually limit the outer radius of a circumstellar disc (or the inner radius of a circumbinary disc). Except 
in the case of a very low-mass companion, the Lindblad resonances are so strong that they are excluded from the disc. In a 
binary star with a mass ratio q = M 2 /Mi of order unity, a circumstellar disc is tidally truncated well inside the innermost 
(m = 2) Lindblad resonance, which in any case lies outside the Roche lobe (Paczynski 1977; Papaloizou & Pringle 1977). For 
smaller mass ratios typical of giant planets orbiting stars, the Lindblad resonances play a more direct role in truncating the 
disc (Lin & Papaloizou 1986). 

Another horizontal tidal effect is the eccentric interaction between the companion and the disc. At eccentric Lindblad 
resonances such as the 3:1 resonance, a local eccentric instability occurs, which may be able to sustain a large-scale eccentric 
distortion of the disc (Lubow 1991). The eccentric corotational resonances act to damp eccentricity, however. 

Studies of the vertical component of the tidal force have mostly been concerned with the case of a companion with 
an inclined orbit. In that case, the vertical tidal acceleration is nearly independent of 2 within the disc. As a result, non- 
axisymmetric bending waves are launched at vertical resonances where the forcing frequency resonates with the natural vertical 
oscillations of the disc (Shu, Cuzzi & Lissauer 1983). In a Keplerian disc the Lindblad and vertical resonances coincide. 

For a companion with a coplanar orbit, however, the vertical tidal acceleration is proportional to 2 and vanishes on the 
mid-plane. Other than those concerned with bending waves, most analytical and numerical studies have relied on a two- 
dimensional treatment of the disc, and the vertical component of the tidal force has been widely neglected. However, Lubow 
(1981) analysed the response of a vertically isothermal disc to vertical tidal forcing by a companion with a circular, coplanar 
orbit. The general response consists of a non-axisymmetric vertical expansion and contraction of the disc, characterized by a 
vertical velocity proportional to 2 . This ‘breathing’ motion corresponds to a compressive mode of the disc, the n = 1 p mode 
in the notation of Lubow & Pringle (1993), or the p) mode in the classification of Ogilvie (1998). At particular radii where 
the forcing frequency resonates with the natural frequency of this mode, a non-axisymmetric p-mode wave may be launched 
and a resonant torque exerted. These are also vertical resonances , but are distinct from those associated with bending waves 
because a different mode is involved. Lubow’s vertical resonances are further from the companion than the corresponding 
Lindblad resonances and, in a thin disc, the resonant torques are much weaker. However, Lubow (1981) pointed out that, 
in a binary star with a mass ratio of order unity, the m — 2 inner vertical resonance typically lies within the standard tidal 
truncation radius, and the associated torque can compete with the viscous torque in the disc. 

Another consequence of the vertical tidal force is a local tilt instability, closely related to the eccentric instability (Lubow 
1992). The tilt instability occurs at inclination resonances that are essentially coincident with the eccentric Lindblad reso¬ 
nances, but the tilt instability is much weaker. 


1.2 Aims of the paper 

In this paper I develop a non-linear theory of Lubow’s vertical resonances. I generalize the analysis of Lubow (1981) to non- 
isothermal discs allowing for non-linearity, viscosity and radiative damping. By reducing the problem to a universal differential 
equation and studying its solutions, I aim to demonstrate the effects of non-linear saturation and dissipation on the resonant 
response of differentially rotating discs to periodic forcing in general. 

In a related paper (Ogilvie 2001, hereafter Paper I) I have presented an analysis of the non-linear tidal distortion of 
a thin, three-dimensional accretion disc by a binary companion on a circular orbit. I showed that the resulting distortion 
can plausibly explain the two-armed non-axisymmetric features seen in the Doppler tomograms of IP Peg and other dwarf 
novae in outburst (e.g. Steeghs 2001). Although the m = 2 inner vertical resonance plays an important role in regulating 
the amplitude and phase of the tidal distortion, the method of analysis in Paper I does not allow for the emission of a wave 
from the vertical resonance. One of the aims of the present paper is to justify that limitation by showing that, under typical 
conditions, vertical resonances are broadened and damped by non-linearity and dissipation to the extent that the emission of 
a wave may reasonably be neglected. 

Before embarking on the detailed analysis, it may be helpful to look ahead to the equation to be derived, 

- +xy+ \y\ 2 y + iby = a. (5) 
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In this equation, x represents the radial distance from the resonant orbit and y the complex amplitude of the response, both 
in dimensionless terms. The linear terms —y" + xy represent the ability of the disc to support freely propagating waves in the 
region x < 0. In the present case these waves correspond to the compressive pj mode. The point x = 0 is the location of the 
resonance, which is also the turning point of free waves. The constant term a on the right-hand side represents the tidal forcing 
in the neighbourhood of the resonance. The inhomogeneous Airy equation — y" + xy = a is familiar in studies of resonant 
wave excitation in differentially rotating discs and features in the analysis of Lubow (1981). The new, non-linear term \y\ 2 y 
derives from couplings between the resonant mode and non-resonant modes. It leads to new effects such as the non-linear 
broadening and saturation of the resonance. The correct derivation of this term accounts for much of the complexity of the 
analysis in the present paper. The linear term i by represents the dissipation of vertical motions by shear and bulk viscosity 
and by radiative damping. It leads to attenuation of the waves and broadening of the resonance. 

The remainder of this paper is organized as follows. In Section 2 I recall the coordinate system and the basic equations 
used in Paper I. In Section 3 I expand the equations in the neighbourhood of the vertical resonance and solve them in a 
systematic manner to derive the complex non-linear Airy equation. In Section 4 I discuss the properties of the equation and 
present numerical solutions for a range of parameter values. The astrophysical consequences of the analysis are explored in 
Section 5. 

2 ANALYSIS 

2.1 Orbital coordinates 

Consider a binary system with two stars in a circular orbit. The stars are approximated as spherical masses Mi and M 2 . Let 
d be the binary separation and 

G(Mi + M 2 ) 
d 3 

the binary orbital frequency. Consider the binary frame rotating about the centre of mass with angular velocity ui, and let 
(R , <j), z) be cylindrical polar coordinates such that star 1 (about which the disc orbits) is at the origin, while star 2 has fixed 
coordinates (d, 0,0). 

I 11 Paper I a system of non-orthogonal orbital coordinates (A, (j> , z) was introduced instead of cylindrical polar coordinates. 
These are based on Paczynski’s orbits for the restricted three-body problem, and naturally account for the principal tidal 
distortion of the disc. The coordinate A is a quasi-radial coordinate that labels the orbits and differs only slightly from R, 
reducing exactly to R in the limit M 2 —> 0. In this paper I again consider the case of a prograde circumstellar disc. However, 
the generalization to retrograde and circumbinary discs is straightforward. 

The standard notation of tensor calculus is adopted, with g ab being the metric tensor, J = yj det(g 0 6) the Jacobian 
of the coordinate system and r/ abc the Levi-Civita alternating tensor. Expressions for the metric coefficients and connection 
components for the orbital coordinate system may be found in Paper I, Section 2. 

2.2 Basic equations 

The basic equations governing a fluid disc in three dimensions are the same as in Paper I. However, the effect of a non-zero 
relaxation time of the turbulent stress will not be considered here, as it was found to introduce algebraic complications without 
affecting the final result of this paper significantly. Accordingly, a purely viscous model of the turbulent stress is adopted here. 


The equation of mass conservation is 

(dt+U a V a )p= -pVaU, (7) 

where p is the density and u a the velocity relative to the rotating frame. The equation of motion is 

p [(ft + u V b )u a + 2rj abc u b u c ] = -pV a (<h + $ t ) - V a p + V b T ab . (8) 

Here oj a is the angular velocity of the rotating frame (only u> z = ui being non-zero), p is the pressure, 

$(77, z) = —GMi(R 2 + z 2 )~ 1/2 - \u 2 R 2 (9) 

is the part of the effective potential not due to star 2, and 

(R,<t>,z) = —GM 2 (R 2 + d 2 - 2dRcos<t> + z 2 )~ 1/2 + GM 2 d~ 2 R cos <j> (10) 

is the tidal potential due to star 2. The turbulent stress tensor T ab is assumed to be given by 

T ab = M ( W + VV) + (^ - | p){X7 c u c )g ab , (11) 
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where p is the effective shear viscosity and pb the effective bulk viscosity. The energy equation is 


( 12 ) 


7 - 1 


(dt + u a S7a)p = - 


7 - 1 


p\7 a u a +T ab V a u b - V a F a , 


where 7 is the adiabatic exponent and F a the radiative energy flux, given in the Rosseland approximation for an optically 
thick medium by 


pa _ 16(77 pjap 

Snp 


(13) 


where a is the Stefan-Boltzmann constant, T the temperature and k the opacity. The equation of state of an ideal gas, 
kpT 

P ~ -j 

f-LuiTTla 

is adopted, where k is Boltzmann’s constant, p m the mean molecular weight and toh the mass of the hydrogen atom. The 
opacity is assumed to be of the power-law form 


(14) 


k = C K p x T y , (15) 

where C K is a constant. This includes the important cases of Thomson scattering opacity (x = y = 0) and Kramers opacity 

(* = 1 , 2 / = - 7 / 2 ). 

The effective viscosity coefficients are assumed to be given by an alpha parametrization. The precise form of the alpha 
prescription relevant to a tidally distorted disc is to some extent debatable. It will be convenient to adopt the form 



where a and ab are the dimensionless shear and bulk viscosity parameters. In the limit of a circular disc in the absence of 
star 2 , this prescription reduces to the usual one, p = up/Cl, etc., where Cl is the orbital angular velocity in the inertial frame. 
For simplicity, the five dimensionless parameters of the disc, (a,ah,'y,x,y) will be assumed to be constant throughout this 
paper. The other parameter in the problem is the mass ratio q = M 2 /Mi. 


3 EXPANSIONS 
3.1 Thin disc 

Analytical progress can be made only if the disc is thin so that certain approximations may be used. The systematic way of 
introducing these approximations is through an asymptotic analysis. 

Let a system of units be adopted such that the radius of the disc and the orbital frequency are 0(1). Let the small 
parameter e be a characteristic value of the angular semi-thickness H/R of the disc. To resolve the vertical structure of the 
disc, introduce a stretched vertical coordinate / according to 

* = eC, (17) 

so that ( = 0 ( 1 ) within the disc. 

Any slow time-dependence of the unperturbed disc, which occurs on a characteristic viscous time-scale a - 1 e - 2 12 -1 , is 
consistently neglected. The desired solution is stationary in the binary frame. 


3.2 Orbital motion 

Let 

= $(R,o) + MR,<t>,o) (is) 

denote the total potential in the binary plane. The equations governing the shape R{ A, 0) of Paczyriski’s orbits and their 
orbital angular velocity f2(A ,</>) can be expressed in the form (Paper I, Section 5.1) 

(R 2 + 2 R% - RR^)n 2 + 2(R 2 + R^uQ = R ^ , (19) 

Qd^[R 2 (Q + ca)] =- , (20) 

where R$ = (dR/d(p)\, etc. These also correspond to the horizontal components of the equation of motion of the disc at 
leading order in e. The divergence of the orbital motion is 
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A = -0 0 (Jfi). 

The primary potential can be expanded in a Taylor series about the mid-plane, 
$ = V 0) (R) + ± e 2 C 2 $ (2) (i?) + 0(e 4 ). 


( 21 ) 

( 22 ) 


The tidal potential can be expanded similarly, but is assumed formally to be smaller by a factor O(e) than the primary 
potential, i.e. 


<f>t = e 


^ (R, 4) + h\ 2 *[ 3) 1 (R, 4>) + 0{V) 


(23) 


Thus <f> = + e^ 1 ^, and Paczyriski’s orbits can be considered as deformed circular orbits such that the orbital quantities 

possess expansions of the form 


R = A + eR (1) (A,</)) + 0(e 2 ), 

(24) 

fi = fi (0) (A) + efi (1) (A, </>) + 0(e 2 ), 

(25) 

A = eA (1) (A,<(>) + 0(e 2 ), 

(26) 

where 



(27) 

Let 


fi = n (0) +w =( G " 1 ) 1/a 

(28) 


denote the angular velocity in the inertial frame. 

An explicit calculation of R^\ and A*- 1 - 1 is deferred until Section 3.8. At present it is sufficient to note that they 
possess Fourier expansions of the form 


_R (1) (A, </>) = Re^ R (1)(m) (A)e iW \ f2 (1) (A, ^>) = Re ^ fi (1)(m) (A)e im *, A (1) (A, 0) = Re ^ A (1)(m) (A)e i 


i m<p 


(29) 


where R (1)(m) and fit 1 )' 771 ) are real owing to the symmetry of the orbits, while A^*- 771 ) is imaginary. 

The assumption that the tidal potential scales as O(e) is a formal device that will turn out to provide a critical level of 
non-linearity in the response. Similarly, it is expedient to choose viscosity parameters that scale formally as 0(e 2 ^ 3 ), i.e. 


2/3 - 
a = t a , 


«b = e 2/3 db- 


(30) 


This will turn out to provide a critical level of dissipation in the response. 

3.3 Expansion about the resonant orbit 

Let A = A* be the resonant orbit, to be identified subsequently, and let a subscript * generally denote a quantity evaluated 
on the resonant orbit. It is known from Lubow (1981) that the characteristic radial extent of the resonance is (H 2 R) 1 ^ 3 , as 
for a Lindblad resonance. One therefore introduces a scaled quasi-radial coordinate £ according to 

A = A, + e 2/3 £, (31) 

so that £ = 0(1) in the resonant region. Any function of A may then be expanded in a Taylor series about the resonant orbit. 
Thus 


R = A, + t 2/3 i + eRi (</>) + 0(e 5/3 ), 

= fl* + e ^ ^ 2 / 3 ^ T efli(<(>) T 0(e ,/3 ), 

A = eAr(<A) + 0(e 5/3 ), 
where 

Mir- 

0 3 /GMiN 1 / 2 

^ 2/3 - ~2a: ’ 

while Ri{(j>) = R^{ A*,0), S2i(<(>) = n^(A*,</>) and Ai(0) = A^(A*,</>). Similarly, 


(32) 

(33) 

(34) 

(35) 

(36) 
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$ (2) = $(2) + e 2 / 3 $®| + e$i 2) (0) + 0(e 4/3 ), 
etc., where = fl* and "hj/i = 2 f)„fl 2 / 3 - 

3.4 Fluid variables 

A steady solution of the basic equations is sought. After careful analysis, the following expansion scheme for the fluid variables 


in the resonant region has been found to be appropriate and self-consistent. 

u x = e 5/3 4 /3 (£,0,C) + O(e 2 ), (38) 

= e 5/3 < /3 (^,0 + O(e 2 ), (39) 

u = e ^ 3 u 4 / 3 (^, </>, C) + e 5 ^ 3 M 5 / 3 (^, <)>, C) + e u 2 (^, <?i, C) + 0(e 7 ^ 3 ), (40) 

P = £ [p*(0 + el/3 Pl/3(£> </>, 0 + e “' /3 / 3 2/3(Ci </>, 0 + e Pl(£i </>! C) + 0( e4,/3 )] > (41) 

P = e s+2 [p*(C) +e 1/3 Pi/3K,<(>,C)+ e 2/3 P2/3(C,<( > ,C) + ePi(?,<(>,0 + O(e 4/3 )] , (42) 

M = e s+8/3 [MC) + e 1/3 Mi/3(^,0,C) + 0(e 2/3 )] , (43) 

Mb = e s+8/3 [Mb*(0 + O(e 1/3 )] , (44) 

T = e 2 [T t (C) + e 1/3 Ti /3 (C,0,C) + O(e 2/3 )] , (45) 

F z = e s+ 11/3 [ F ; (C ) +£ l/3^ /3(? , 0;C ) + O ( e 2/3 ) j (46) 


Here s is a positive parameter, which drops out of the analysis, although formally one requires s = (10 — 6 y)/(6 + 3x) in order 
to balance powers of e in the opacity law. Other components of F a are smaller and will not be required. 

According to this scheme, the relative amplitude of the perturbation in the resonant region is O^e 1 ^ 3 ). This is the 
magnitude of the fractional perturbations of density and pressure, and of the Mach number of the vertical velocity perturbation. 
The relative amplitude of the perturbation is greater than the relative magnitude of the tidal potential, O(e), because of the 
effect of the resonance. The horizontal velocity perturbations are smaller than the vertical one by a factor O^e 1 ^ 3 ). 

3.5 Equations to be solved 

The above expansions for the fluid variables are now substituted into the basic equations, which are expanded only as far as 
is required for the analysis that follows. 

3.5.1 Equation of mass conservation 
0 {e s+1/3 ): 

^9 ( />Pi/3 + u Z 4/ 3 d<; P* = ~p*d c u i/3 . (47) 

0(e s+2/3 ): 

^l*d ( j > p2/3 + w|/ 3 dfPi / 3 + u 3 / 3 9^p* = —p*d^U 5/ 3 — pi/ 3 d(ul/ 3 . (48) 

0(e s+1 ): 

Q.td^pi + «4/ 3 dcP2/3 + (^2/3^9,/, + U5/ 3 df)pi/ 3 + u\dcp* = —p*(d^u 3 / 3 + + Ai) — Pi/ 3 9quI/ 3 — P2/3^c' u 4/3- (49) 

3.5.2 Energy equation 
0(e s+7/3 ): 

0*d^pi /3 + ul /3 d ( p* = - 7 p*< 9 C M 4 /3 . (50) 

0(e s+8/3 ): 

^*9^2/3 + u z 4/3 9 c p 1/3 + ul /3 9t;p, = - 7 p*d C Mg /3 - 7 Pi/ 3 d<;uI /3 + (7 - 1) \ /3 - 9 c Ff) . (51) 

0(e s+3 ): 

^1*9^1 + ul/ 3 9 C p 2 /3 + (^ 2 / 3 ?d 0 + ul /3 9 C )p 1/3 + U Z 29 q p* = - 7 P*(<%M 5 /3 + 9 c u\ + Ai) - 7 Pl/ 3 ' 9 C « 5/3 - 7 P 2 / 3 ^c M 4/3 

+(7 - 1) (pi/ 3 A 2 ^ /3 - d ( Ft /3 ) . (52) 
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3.5.3 Equation of motion 
A-component, 0(e s+5 ^ 3 ): 
p*{£l*d(j)U3/3 2A*r2*iig/g) — cl^Pi/3- 
(/(-component, 0(e s+5/ ' 3 ): 

P* ^ m 5/3^2/3 + ^»^ M 5/3 + = 

2 -component, 0(e s+1 ): 

0 = —p»<&S 2) C “ dqP*- 
2 -component, 0(e s+4//3 ): 
p*n*a0M4 /3 = -Pi/3$i 2) C - <9fPi/3- 
2 -component, 0(e s+5 / 3 ): 

P* (fi,a^M5 /3 + M4/3^C«4/3) + pl/3^*d<t>Ul /3 = -p 2 /3$? ) ( - p.$2/UC ~ ^CP2/3- 
2 -component, 0(e s+2 ): 

P* [f2*90U2 + M4/3^M5/3 + (f2 2 /3^5^ + W5/3<?c) M 4/3] + Pl/3 (f2*9<j>U5/ 3 + W4/3^C u 4/3) + P2/3^*^<A u 4/3 

= -pi$i 2) C - Pl/3*2/UC - P*$i 2) C - P*^C - d( [pi - (pb* + |p*)9cM4/3] ■ 


3.5.4 Constitutive equations 
16aT*~ y 


Ff = — 


P* = 


3C K p 
k p^'T* 


l-fx 




ap* 



pb* — 


QbP* 

n* 


Also 

^ = {3 - y)I rr 


Pl/3 

p* 


Pl/3 t 1/3 

p, t* ’ 


(1+ *)*£ + 


dCCi/3 

d c T* ’ 


&P1/3 

111,3 - ir~ 


3.6 Solution 

3. ff. i Unperturbed vertical structure 

From equation (^5|), one has the important relation 

d c p* = -pMlC 

which expresses the vertical hydrostatic equilibrium of the unperturbed disc. The final two terms in equation ( |bI| ) 
taken to balance, since these represent the unperturbed thermal equilibrium of the disc, i.e. 

d c F: = | apM*. 

3.6.2 Form of the resonant mode 

Equations (^), (^) and (^) may be combined to give 

Cut/ 3 = 0, 

where 


(53) 

(54) 

(55) 

(56) 

(57) 

(58) 

(59) 

(60) 

( 61 ) 

(62) 

(63) 

(64) 

(65) 
may be 

( 66 ) 

( 67 ) 
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£ = p.(#a2+n;)-0 f (7P.0 c ) (68) 

is a linear operator. Normally £ would be a non-singular operator and equation (^57|) would have no non-trivial solution. 
However, the resonant orbit is, by definition, one on which £ is a singular operator possessing a null eigenfunction of the form 


17710 > 

w = e 

such that Cw = 0. This requires 
m 2 fl 2 t = (7 + 


(69) 


(70) 


which is the condition given by Lubow (1981) for a vertical resonance for azimuthal wavenumber m (m > 0 always). Accord¬ 
ingly, the general solution of equation ([si]) is 

u 4/3 = R- e e lm ^ imfi./C] (71) 

where /(£) is a dimensionless function to be determined. The corresponding solutions of equations ( ffTj ) and (|H(]) are 
Pi/s = Re [e lm ^ /(p« + C^c P*)] > (72) 

Pi/3 = Re [e' m ^ /(7P* + C^cP*)] • (73) 

One also finds 

Ti/a = Re {e im< ^ / [(7 - 1 )T. + C^T.]} , (74) 

* 1/3 = Re {e im0 / [((4 - y)( 7 - 1) - *)*? + C^]} . (75) 

This motion represents a simple ‘breathing mode’ of the disc. In the notation of Ogilvie (1998) it is the pj mode, i.e. the first 
p mode of even symmetry about the mid-plane. 


3.6.3 Solvability condition for the linear operator 
At higher orders one naturally obtains equations of the form 

CX = F, (76) 


where X is unknown and F is known. Now £ is a singular operator, and it is important to understand the conditions under 
which such an equation is soluble. Since £ is self-adjoint, the solvability condition is easily obtained in the form 


w*F d(/>dC, = 0, 

where w is the null eigenfunction (|69|), and the integration is over the full azimuthal and vertical extent of the disc. 


(77) 


3.6.4 Horizontal velocities 

Equations (^) and ( |54| ) may be combined to give 
p*(ft*d| + &*)u 3/3 = 


The solution is 


^5/3 


= Re 


im(j} imfl* df 1 , -r, . 

e v —— — — (7P* +C%M 


7^* d £ P* 


m 5/3 = Re 


i m<t> 


d / 1 


27A* fl* d£ p 

where equation (^) has been used. 


(7P* + P*) 


3.6.5 Excitation of the second harmonic 
Equations and (0 may be combined to give 

•£^5/3 = F 5/3, 

where 

F 5/3 = —^l*d^(p*ul/ 3 d(ul/ 3 + Pi/3fl*9^M4/ 3 ) + &l({ul/ 3 d(Pi/ 3 + Pl/3<9c M 4/ 3 ) + ^c( u 4/3^CPl/3 + 7Pl/3^C M 4/3 ) • 


(78) 


(79) 

(80) 


(81) 

(82) 
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This evaluates to 


F 5/3 = Re 


2i m<t> 


;(7 + 1)(7 + 3)i 


The solvability condition is satisfied, and the solution is 


u 5/3 — Re 


2i m<f) 


6 


(7 + 3)i mQtjC 


The corresponding solutions of equations (Etq) and (pTI) are 


P 2/3 = Re je 
P 2/3 = Re {e 2im?1 i/ 2 


2i rntfi 


(7 + 3)(p* + C^cp*) + (p* + C^cp*) + C<9c(p* + C^cP*)J | + P2/3 (^C)j 
(7 + 3)(7P* +C^CP*) +7(7P* + C^CP*) + C9c(7P* +C^CP*)] | + P^/ 3 (£,C), 


where the axisymmetric parts satisfy 

0 = -P2/3$i 2) C - P*<Uc - «3 - i(7 + i)^|/| 2 a c (p,c 2 ), 


(83) 

(84) 

(85) 

( 86 ) 

(87) 


and represent the Taylor expansion of the unperturbed vertical structure about the resonant orbit, plus a non-linear term due 
to the wave. 


3.6.6 Excitation of the resonant mode 

Equations ([fi]), (^) and (^i|) may be combined to give 

Cu z 2 = F 2 , 

where 

F 2 = — f2*<9<£ jp* [u4/ 3 9£lig/3 + (n 2 /3 ^dcj, + Ug/39c) w 4/3] + Pl/3 (H*^Ug/ 3 -f- W4/3^C U 4/3) + P2/3^*^0 M 4/3 

+Pl/3^y 3 ^C + P*$ 1 2) C + P^i^C ~ d ( [(pb* + |p*)9cW4/ 3 ] I 

[«4/39cP2/3 + (^2/3 + «5/3<9c)Pl/3 + P*(<9£%/3 + ^l) + Pl/3^C M 5/3 + P2/3'9c«4/ 3 ] 

+di [«4/3%>2/3 + {^ 2 / 3 ^ + W5/3^c)Pl/3 + 7P*(<%«5/3 + A l) + 7Pl/3dc w 5/3 + 7P2/3<9c M 4/3 

-(7-1) (pi/ 3 A 2 ,^ /3 - Sfitf/s)] ■ 

This equation need not be solved in detail. It is sufficient to apply the solvability condition, 



c iTa d^dC = o. 


which yields the required equation for /(£). After multiplication by i/vr, this has the form 

Ci0 - C 2 U - i C 3 f - C 4 |/| 2 / + Cs = 0, 

where the five real, constant coefficients are given by 

Cl = -7ZT [— (7P*+ C9cP*) 2 dC, 

7 H 2 J p * 

c 2 = 2(7 + i)mn„(n* — f)»)n 2 /3 / p*C 2 dC, 


c 3 = | (7+ l)(«b + fa) + | (7 - 1) [(3 - y)(7 - 1) - *]a} J P*£ 2 dC, 
C 4 = ^_7(7 + l)(3-7)m!7»fi2 J P*C 2 dC, 


C5 = mil, 




(2)(m) 


+ $ 


(3)(m) 


J p*C 2 dC+( 7 -1) ii* J P*C 2 dC- 


Equation Q91[) may be reduced to a normal form, 


- + \y\ 2 y + % = a, 


( 88 ) 


(89) 

(90) 

(91) 

(92) 

(93) 

(94) 

(95) 

(96) 

(97) 
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10 G. I. Ogilvie 


by means of the transformations 


x = c-r 1/3 cn, 

(98) 

y( x) = c; 1/6 c~ 1/3 cl /2 f(0- 

(99) 

The two essential parameters are 


a = C'" 1 / 2 C 2 - 1 C 4 1 / 2 C , 5 , 

( 100 ) 

b = C~ 1/3 C~ 2/3 C 3 , 

( 101 ) 

and are both real. 



3.7 Relation to physical units 


When the asymptotic scalings are removed, the dimensionless parameters can be expressed directly in terms of the physical 
variables evaluated at the location of the resonance. The ordering parameter t then cancels out. However, the integral involved 
in Ci introduces the dimensionless parameter e defined by 

2 ___./■ _ 

( 102 ) 


1 / . dp\“ , . 2 , 2^4 f 2 , 

- yyp + z J dz = e A*12, / pz dz. 


This is a more specific measure of H/R, and depends on the vertical structure of the disc. For a vertically isothermal disc, 
e 2 = (y 2 — 27 + 3 )(H/R) 2 , where H is the isothermal scale height. For a vertically polytropic disc with polytropic index 3/2, 
e 2 = (l/40)(77 2 — IO 7 + 15 )(H/R) 2 , where H is the true semi-thickness. For a vertically polytropic disc with polytropic index 
3, e 2 = (1/44)(5 7 2 - 87 + 12 )(H/R) 2 . 

A quantity deriving from C 2 is the rate of detuning of the vertical resonance, 
d 


D = 


din A 


(7 + 1)S7 2 — m 2 n 2 . 


Another quantity having the dimensions of frequency-squared is the forcing combination appearing in C 5 , 


4- = )(m) + $® (m) + (7 - 1 ) 


A 


(m) 


imfT 


n 


One then finds 

7 2 (7+ 1)(3 — 7 ) 


12 


1/2 


JL 

15' 


This shows why the tidal potential was chosen to scale formally as O(e), in order to have a = 0(1). Similarly, 


b = |(7 + 1 )( a b + fa) + |(7 - !) [(3 - y )(7 - 1) - *] a| 


'yQt 


1 1/3 


(7 + l) 3 / 2 e 2 X > 2 


The scaling between x and the physical length A is given by 
A,. 


dA 

dx 


i 2 nl ^ 1/3 


7 ® 


(103) 


(104) 


(105) 


(106) 


(107) 


3.8 Evaluation of the forcing terms 

There are three contributions to the forcing coefficient 4/. The direct vertical forcing involves two terms: which is the 

perturbation of the primary potential due to the tidal distortion of the orbit, and which is the tidal potential itself. 

The term proportional to is an indirect forcing of the p mode resulting from the divergence of the tidally distorted 

orbital motion. 

In the following it is assumed that the resonance is an inner vertical resonance with m > 2. The Fourier components of 
the tidal potential can be expressed in terms of the Laplace coefficients 


6^ m ^(/3) = — / (1 + (3 2 — 2/3cos <f>) 7 cos(m^) d<j>, 


(108) 


where (3 = \/d. Thus 

^ )(m ) = &$(/?), <f>f )(m) = b( $W- 


(109) 


The location of the resonance is 
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Vertical resonances in accretion discs 11 


^ V3 h - -L (7 + 1) 1 /2j 2/3 _ 

a m J 

The orbital equations (^9[) and (^) yield the following linear equations for the perturbed orbital quantities: 
(m 2 fl (0)2 + SO. 2 ) AR (1)(m) + 2fiA 2 fl (1)(m) = , 

A12 (0) (Af2 (1)(m) + 21LR (1)(m) ) = -$[ 1)(m) . 

Thus 

(m 2 n {0)2 - Cl 2 ) A7? (1)(m) = A^$l 1)(m) + 

and 

A (1)(m) = i m [fi (1)(m) + n (0) ( fl(1)(m) + — R (1)(m A 

A dA J 


1 ^(l)(m) 2 Cl (1 )( m ) o(0 ) / 

'Him 1 X K ' 


^ (1)(m) , indlW 

A dA 


* 


Evaluating these quantities at the resonance, one finds 

(2 Km) _ g GMl p(m) 


Ai 




,(m) _ f ^(SRm) 30* 

t* ' \ 2 q 

A* a 


and therefore 
GM\ q 


® = 


d 3 7 # 


(Z 3 *) + 3m (7 + !) 1/2 (7 - l)^*i/2 (A 


) 


+ T [ 2 7 - 1 - m(7 + !) 1/2 (7 - 1)] P*b ( ™o'(P*) + 2m(7 + 1) 1/2 &i 72 (Z 3 *) 


Here the condition for an inner vertical resonance, 
m. 0 , = (7 + l) 1 / 2 /), 

has been used, as well as the fact that the tidal potential satisfies Laplace’s equation, 

1 _d_ 

AdA 

One also Hnds 

GMi 


^A-j^X" 0 ) _ $«(”•) + $(3)(m) = o 


V = 


d 3 


[3m( 7 + l) 1/2 (l + 9 ) 1/2 X 3/2 ]. 


( 110 ) 

( 111 ) 

( 112 ) 

(113) 

(114) 

(115) 

(116) 
(117) 


(118) 


(119) 


( 120 ) 


( 121 ) 


3.9 Illustrative parameters 


The following parameter values might be appropriate for a system such as IP Peg: 7 = 5/3, m = 2 and q = 0.5. For the 
purposes of illustration, suppose that Ob = 0, x = 1 and y = —7/2. Then one Hnds 

A* « 0.282 d, a « 0.0110 e -1 , b « 4.20 ne _2/3 , ^ « 0.196 e 2 / 3 d. (122) 

ax 

For these parameters, radiative damping provides a somewhat greater contribution to b than does viscous damping. A more 
non-linear and less dissipative example occurs when 7 = 6/5, m — 2 and q = 1: 

A. « 0.322 d, a « 0.0725 e -1 , b « 1.26 ne _2/3 , ^ « 0.230 e 2 / 3 d. (123) 


4 ANALYSIS OF THE COMPLEX NON-LINEAR AIRY EQUATION 
4.1 Physical interpretation 

The non-linear ordinary differential equation ([)?]) describes the response of the disc in the neighbourhood of the vertical 
resonance. All the fluid variables can be determined from the solution of this equation. A detailed study of the solutions is 
therefore called for. 
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Figure 1. Mode couplings producing the cubic non-linearity. Quadratic non-linear terms in the basic equations, such as u z d z u z , couple 
together the different azimuthal mode numbers. The combination of any two modes mi and m 2 (shown by arrows converging on a dot) 
forces modes |mi ±m 2 | (each shown by an arrow leaving a dot). Although all modes are forced at similar levels by the tidal potential, the 
resonant mode m responds with anomalously large amplitude. Its non-linear self-interaction excites modes 0 and 2 m at larger amplitudes 
than they would acquire by direct forcing, and their own interactions with mode m feed back on mode m. The net result is that mode 
m affects itself through a cubic non-linearity. 


Equation ( |97| ) does not appear to have been studied before. It may be conveniently described as the complex non-linear 
Airy equation. The five terms may be interpreted as follows. The linear terms —y" + xy represent the ability of the disc to 
support freely propagating waves in the region x < 0. In the present case these waves correspond to the compressive pf mode. 
The point x = 0 is the location of the resonance, which is also the turning point of free waves. The non-linear term \y\ 2 y 
derives from mode couplings. The linear term i by represents the dissipation of vertical motions by shear and bulk viscosity 
and by radiative damping. The term a is the tidal forcing. 

It is important to realize that the cubic form of the non-linear term is not a low-order truncation of a more complicated 
expression, but is the exact form of the non-linearity within the consistent ordering scheme described above. The cubic term 
arises as follows. For simplicity, consider only the vertical velocity u z , which is proportional to z. In the equation for u z there is 
a non-linear term u z d z u z , also proportional to z, that couples together different azimuthal modes. With the scalings adopted 
in this paper, the vertical tidal acceleration —d z &t is 0(e 2 ). All azimuthal modes respond to this forcing with u z = 0(e 2 ), 
except the resonant mode m, which responds with u z = 0(e 4 ^ 3 ). The quadratic self-interaction of mode m forces modes 0 and 
2 m at 0(e 5 / 3 ), i.e. at lower order than the tidal forcing, and they respond with u z = 0(e 5 ^ 3 ). The interaction of these modes 
with mode m forces modes m and 3m at 0(e 2 ), i.e. at the same order as the tidal forcing. This modifies the leading-order 
resonant mode, leading to a cubic non-linearity (Fig. 1). 


4.2 Special limits 


When a <C 1, the response is y = O(a) and the cubic term may be neglected. One then has 


dx 2 


+ xy + i by = a. 


In the absence of dissipation, b = 0 and one obtains an inhomogeneous Airy equation, 


d 2 _y 

dx 2 


+ xy = a, 


which is well known in studies of resonant wave excitation in discs. The general solution is 
y = iraG\{x) + ciAi(a;) + c2Bi(a;), 


(124) 


(125) 

(126) 


where ci and C2 are constants, Ai and Bi are the Airy functions, and Gi is the inhomogeneous Airy function defined by, e.g., 
Abramowitz & Stegun (1965). As x —> +00, both Gi(a;) and Ai(:r) decay monotonically, while Bi(*) diverges exponentially. 
Therefore C2 = 0 for an acceptable solution. For x < 0, both Gi(:r) and Ai(a;) are wavelike. The desired solution is 

y = na [Gi(a:) — iAi(ac)], (127) 


which has the asymptotic form 
y ~ 7r 1 ^ 2 a(—a:) -1 / 4 exp 


2., ,3/2 w 

3 1 ~ x T 


(128) 


as x —* —00, and represents a purely outgoing wave emitted from the resonance and travelling to the left. For any other choice 
of ci, the solution would contain an incoming wave component that is unphysical unless the outgoing wave is able to reach 
an edge of the disc, reflect and return to the site of launching with a non-negligible amplitude. 
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When dissipation is present so that 6 > 0, the emitted wave is damped. Formally, the solution may be obtained simply 
by replacing x by x + i b. This effectively moves the resonance off the real axis. 


4.3 Angular momentum flux and tidal torque 

One can multiply equation ( J97[ ) by y* and take the imaginary part to deduce the conservation law 

^ im (- y *S) =aim( - y) - %|2 - (129) 

In the absence of dissipation and forcing, there is a strictly conserved quantity which is proportional to the angular momentum 
flux of the wave. The term alm(-y) is proportional to the tidal torque density, and the term — b\y\ 2 represents dissipative 
losses. The total tidal torque exerted on the disc is therefore proportional to 


/ OO 

Im(— y) Ax. 

-OO 

In the linear, inviscid limit (a<l, 6 = 0) this quantity is equal to 


/ Ai(a:) da: = iva 2 . 


J — OO 

Define the tidal torque parameter 

1 f°° 

fr = — / Im(—y) da:. 

7ra / 

«/ —OO 


(130) 


(131) 


(132) 


In linear theory f T = 1. When 6 = 0 the torque goes entirely into launching the wave. When dissipation is included, it is still 
true that f T — 1, because the contour integral is unchanged when the function is shifted parallel to the imaginary axis. The 
torque then goes partly into launching the wave and partly directly into changing the angular momentum of the disc. As the 
wave attenuates, its angular momentum is also transferred to the disc. In non-linear theory /r^l in general. 

It can be shown that the tidal torque, in physical terms, is equal to 


T = —f T 


7rA“mT 

V 


pz 2 d z, 


(133) 


where all quantities are to be evaluated at the resonance A = A*. This is closely related to the standard formula for the torque 
exerted at a Lindblad resonance (Goldreich & Tremaine 1979) except that (i) vertical rather than horizontal forcing is involved, 
(ii) the second vertical moment of the density, rather than the surface density, appears, and (iii) the denominator involves 
the rate of detuning of the p-mode resonance condition rather than the epicyclic resonance. The factor fr also accounts for 
a non-linear correction to the torque. In the case fr = 1 and for an isothermal disc with 7 = 1, this formula agrees with 
equation (40) of Lubow (1981). 


4.4 Numerical method of solution 


In the non-linear case, equation (^?j) must be solved numerically. In principle, it might be solved as a two-point boundary-value 
problem with suitable conditions at large positive and large negative x. However, no satisfactory method was found to solve 
the non-linear equation by this method owing to the existence of unwanted growing solutions. 

Instead, equation (|97j) has been solved by converting it into a time-dependent partial differential equation, 


.dy d 2 y 

+xy + 


y + iby = a. 


(134) 


dt dx 2 

The steady solutions of this equation are identical to the solutions of equation (]i)7|), and can be found numerically by solving 
the time-dependent equation as an initial-value problem starting from a suitable initial condition such as y = 0. In fact, this 
is not a purely formal device. It is possible to extend the analysis of this paper to permit the solution to vary slowly in time 
in the binary frame. This can be done by allowing the perturbations to depend on a slow time coordinate e 2 / 3 f in addition to 


£, cj> and £. The resulting amplitude equation is precisely equation (134), with a suitably scaled dimensionless time variable. 


Therefore time-dependent solutions of equation (134) that approach a steady state genuinely describe the development of a 


steady wave p atter n in an initially undisturbed disc, and also suggest the stability of the steady solution. 

Equation (134) is a dispersive, non-linear wave equation related to the non-linear Schrodinger equation and the Ginzburg- 
Landau equation. The term xy breaks the translational invariance, providing a turning point for waves at x = 0. The term 
iby with 6 > 0 ensu res that all solutions decay in the absence of the forcing term a. 

Equation (134) has been solved by discretizing it in space on a regular grid of 5000 points on the interval —100 < x < 100. 
A second-order, centred finite difference approximation of the second derivative was used. The desired solution is evanescent 
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for large positive x and has the form of an outgoing wave at large negative x. In order to select this solution and to prevent 
wave reflection as far as possible, the computational domain was extended to —110 < x < 110 with an additional critical 
damping effect applied in the added intervals. The boundary condition y — 0 was then applied at x = ±110. The resulting 
ordinary differential equations in the temporal domain were solved using a fifth-order Runge-Kutta method with adaptive 
step-size. Starting from the initial condition y — 0, the solution was advanced until a steady solution was obtained. The torque 
parameter was then evaluated by integrating the solution over —100 < x < 100 . 

Some illustrative solutions of the complex non-linear Airy equation are shown in Figs 2 and 3. For each solution, the 
estimated value of the torque parameter fr is quoted. It should be borne in mind that, for those solutions in which a wave of 
significant amplitude reaches the edge of the computational domain at x — — 100 , the integrated torque continues to fluctuate 
spatially at the level of a few per cent, and so the value of fr depends to this extent on exactly where the boundary is drawn. 

In Fig. 2 the effect of increasing the dissipation parameter b on solutions in the linear regime a 1 is investigated. For 
6^1 the most noticeable effect is the attenuation of the emitted wave. For 6 <; 1 the emitted wave is almost completely 
absent and the resonant response is a broad, smooth distortion. The estimated torque parameters are consistent with the 
analytical expectation that f T — 1 when a<l, for any value of b. 

In Fig. 3 the effect of increasing the amplitude parameter a is investigated. The most obvious effect is a broadening of 
the resonance, in which the first wavelength of the emitted wave becomes distorted into a much wider feature. The amplitude 
of the emitted wave increases initially more rapidly than a, and fr achieves values in excess of unity. For larger a, however, 
fr declines to values less than unity, perhaps indicating a saturation of the resonance. 

When either a or b is large, the resonance is substantially broadened and the scaling assumptions used to derive the 
complex non-linear Airy equation become inapplicable. Under these conditions, the vertical resonance can be treated, as in 
Paper I, as a part of the smooth tidal distortion of the disc. 


5 DISCUSSION 


In this paper I have developed a general theory of vertical resonances, first analysed by Lubow (1981), which are an important 
aspect of the interaction between an accretion disc and a massive companion with a coplanar orbit. When Lubow’s analysis 
is generalized to allow for non-linearity of the response, and dissipation by radiative damping and turbulent viscosity, the 
problem reduces to a universal, non-linear ordinary differential equation with two real parameters. 

Numerical solutions of a time-dependent version of the complex non-linear Airy equation describe the process by which a 
steady and stable pattern of non-axisynmietric distortion is established in an initially undisturbed disc. For small values of the 
dissipation parameter b, a p-mode wave is launched at the resonance and propagates radially through the disc, carrying angular 
momentum. For larger values of b the wave is attenuated or not launched at all, and the resonant response is broadened into 
a smooth feature. The effect of increasing the amplitude parameter a is to broaden the resonance and ultimately to saturate 
it. The total tidal torque exerted at the resonance is independent of b in the linear regime aCl, but can be either greater 
or less than the linear prediction when non-linearity is important. 

If the magnitude of the resonant tidal torque exceeds that of the local viscous torque, the disc may be truncated at the 
resonance (Lubow 1981). The viscous torque in the absence of tidal distortions is 


d^dz = — 37raA 2 f2 2 / pz 2 dz, 
qA / 


and a comparison with equation (132) shows that truncation occurs if a < atrunc, where 
/tuiT 2 


C^trunc 


wn 2 


(135) 


(136) 


However, if a is so small that the launched wave does not damp significantly in the vicinity of the resonance, the tidal torque 
is not transmitted to the disc locally and the process of truncation cannot occur in a straightforward manner. In the derivation 
of the complex non-linear Airy equation it was assumed that a = 0(e 2 ^ 3 ) while T = O(e). Therefore the process of tidal 
truncation is not described by that equation, but occurs in a different parameter regime. 

Some illustrative values of atrunc for the m = 2 inner vertical resonance are given in Table 1, where it is assumed that 
/t — 1. It appears that truncation by this resonance is unlikely in cataclysmic variable discs, even in quiescence, and especially 
if 7 « 5/3. Although Lubow (1981) argued that truncation is likely, his estimates were based on the case 7 = 1 , for which the 
resonant radius is maximal. This results in a much larger torque than is obtained when 7 = 5/3, because of the proximity to 
the companion. However, for the same reason, the resonance is less likely to lie inside the disc in the case 7=1. 

Much more probable is that a non-destructive, two-armed structure is formed, as described by the complex non-linear 
Airy equation. In Paper I I have argued that this effect, in the limit in which the resonant response is broadened and damped 
by non-linearity and dissipation, could explain the non-axisymmetric structures seen in Doppler tomograms of dwarf novae in 
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Table 1 . Critical value of the viscosity parameter below which the disc may be truncated at the m = 2 inner vertical resonance. 


9 

7 

^trunc 

1 

1 

0.07489 


6/5 

0.01694 


5/3 

0.00037 

0.5 

1 

0.03533 


6/5 

0.00796 


5/3 

0.00018 

0.2 

1 

0.00935 


6/5 

0.00210 


5/3 

0.00005 


outburst. Indeed, estimates of the parameters for a system such as IP Peg (Section [3.S| ) suggest that the resonance is strongly 
damped and mildly non-linear. 

It is likely that the techniques used in this paper will be useful in other circumstances where spatial resonances occur in 
wave-bearing media subject to periodic forcing. The equation derived has a concise form and its terms have a clear physical 
interpretation, suggesting that it may be generic to some extent. In any case it provides a simple mathematical model in 
which the effects of non-linearity and dissipation on the resonant launching of waves can be studied. 
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